A Method for Measuring Elliptic Flow Fluctuations in PHOBOS 
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We introduce an analysis method to measure elliptic flow (112) fluctuations using the PHOBOS 
detector for Au+Au collisions at ^/sNN = 200 GeV. In this method, V2 is determined event-by- 
event by a maximum likelihood fit. The non-statistical fluctuations are determined by unfolding the 
contribution of statistical fluctuations and detector effects using Monte Carlo simulations(MC). Ap- 
plication of this method to measure dynamical fluctuations embedded in special MC are presented. 
It is shown that the input fluctuations are reconstructed successfully for (1)2) > 0.03. 



I. INTRODUCTION 



Studies of collective flow have proven to be fruitful 
probes of the dynamics of heavy ion collisions. Ellip- 
tic flow (V2) in heavy ion collisions was first discussed in 
[l)] and has been measured at AGS^^, and SPS0, [1| 
energies. The first measurement of elliptic flow at RHIC 
was performed by the STAR collaborationfs'l. PHOBOS 
has measured elliptic flow as a function of pseudorapid- 
ity, centrality, transverse momentum, center-of-mass en- 
ergy and nuclear species 0, H, H, In particular, the 
measurements of V2 as a function of centrality provide in- 
formation on how the azimuthal anisotropy of the initial 
collision region drives the azimuthal anisotropy in par- 
ticle production. When two nuclei collide with non-zero 
impact parameter, the almond-shaped overlap region has 
an azimuthal spatial asymmetry. If the particles do not 
interact after their initial production, the asymmetrical 
shape of the source region will have no impact on the 
azimuthal distribution of detected particles. Therefore, 
observation of azimuthal asymmetry in the outgoing par- 
ticles is direct evidence of interactions between the pro- 



duced particles. In addition, the interactions must have 
occurred at relatively early times, since expansion of the 
source, even if uniform, will gradually erase the magni- 
tude of the spatial asymmetry. Hydrodynamical models 
can be used to calculate a quantitative relationship be- 
tween a specific initial source shape and the distribution 
of emitted particles fll||. At the high RHIC energies, the 
elliptic flow signal at midrapidity in Au-|-Au collisions 
is as large as that calculated under the assumption of 
a boost-invariant relativistic hydrodynamic fluid. The 
presence of a large flow signal has been considered to be 
a proof of early equilibration in the colliding system p^. 

The azimuthal anisotropy of the initial collision region 
is quantified by the eccentricity of the overlap region of 
the colliding nuclei. The customary definition of eccen- 
tricity, which we call the "standard eccentricity," is deter- 
mined by relating the impact parameter of the collision 
in a Glauber model simulation to the eccentricity calcu- 
lated assuming the minor axis of the overlap ellipse to be 
along the impact parameter vector. Thus, if the x axis is 
defined to be along the impact parameter vector and the 
y axis perpendicular to that in the transverse plane, the 
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eccentricity is defined by: 



^standard 



(1) 



where ax and (Jy are the RMS widths of the participant 
nucleon distributions projected on the x and y axes, re- 
spectively. 

Measurement of V2 fluctuations as a probe of early 
stage dynamics of heavy-ion collisions has been sug- 
gested earlier by Mrowczynski and ShuryakpLSj- How- 
ever, Miller and Snellings have pointed out that fluctua- 
tions in the shape of the initial collision region must be 
understood first before addressing other physical sources 
of V2 fluctuations[l3|. Furthermore, the latter argue that 
experimental measurements of V2 can be affected by these 
fluctuations. They show approximate agreement between 
the predictions from a fluctuating eccentricity model and 
the differences in V2 measures obtained via two, four 
and six particle cumulant methods in Au-f Au collisions, 
where the standard definition of eccentricity is used. 

The elliptic flow in Cu-f-Cu collisions is observed to 
be surprisingly large, particularly for the most central 
eventsflO]. PHOBOS has proposed that event-by-event 
fluctuations in the shape of the initial collision region 
can be a possible explanation for the large V2 signal in 
the small Cu-|-Cu system[lol]. For small systems or small 
transverse overlap regions, fluctuations in the nucleon po- 
sitions frequently create a situation where the minor axis 
of the ellipse in the transverse plane formed by the par- 
ticipating nucleons is not along the impact parameter 
vector. An alternative definition of eccentricity, called 
the "participant eccentricity", tpart, is introduced to ac- 
count for the nucleon position fluctuations such that the 
eccentricity is calculated with respect to the minor axis 
of the ellipse defined by the distribution of participants 
found using a Monte Carlo approach. Using the same 
coordinates as before: 



4cr2 



^part 



(2) 



where cr-r 



'xy — {xy) — {x){y). The average values of 
^standard and Spart arc quitc similar for all but the most 
peripheral interactions for the Au+Au system. For the 
smaller Cu-|-Cu system, however, fluctuations in the nu- 
cleon positions become quite important for all central- 
ities and the average eccentricity can vary significantly 
depending on how it is calculated [lOj. 

Understanding fluctuations in the initial collision re- 
gion has proven to be crucial to interpret the (U2) re- 
sults. Measurement of V2 fluctuations can be used to test 
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FIG. 1: Participant eccentricity, epart, of the collision zone 
as a function of number of partipicant nucleons, Npart, for 
Au+Au collisions at ^/sWn = 200 GeV from the PHOBOS 
Glauber MC. The black points show the average epart O. 
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FIG. 2: Squares show the average elliptic flow parameter, V2, 
measured near mid-rapidity, as a function of Npart for Au+Au 
collisions at ^/sJfN = 200 GeV. Only statistical errors are 
shown. Dots show the RMS of event-by-event V2 distribution 
from PHOBOS Glauber MC (see text.) 



the participant eccentricity model. Fig [T] shows a distri- 
bution of Glauber model simulated events as a function 
of Spart and centralityjlOj. These simulations have been 
used to caculate the mean, i^^part) and the RMS, depart ^ 
of participant eccentricity as a function centrality. Ideal 
hydrodynamics leads to V2 oc efll]. Assuming this holds 
event-by-event, this condition would imply that: 



(3) 



crv2 _ 
V2 e ' 

where ay^ is the RMS of the event-by-event distribution 
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of V2. PHOBOS (■U2) results (Fig. ^ can be used to 
estimate cr^^ by Eq. [3l These estimates are also shown 
in Fig.[2]for Au+Au collisions at = 200 GeV. It is 

important to note that in these estimates, other sources 
of elliptic flow fluctuations are neglected. 

II. METHOD 
A. Overview 

We considered three different methods to measure el- 
liptic flow fluctuations: 

• Measuring (ti|) via two particle correlations and 
extracting V2 fluctuations by a comparison to the 
(V2) results. 

• Measuring event-by-event. 

• Measuring V2 event-by-event. 

Two particle correlations in AA collisions can be used 
to measure u| [lB| ■ V2 fluctuations can be calculated com- 
paring (wl) and (U2) by: 

< = («2> - {V2r (4) 

However, the (^2) measurement has significant sys- 
tematic uncertainties, mainly due to the uncertainties 
in the reaction-plane determination [9^. The mea- 
surement will have systematic uncertainties from various 
other sources. The comparison of two quantities obtained 
with different techniques will make it hard to extract 
precise results. A similar approach has been suggested 
in [1^, where the effects of azimuthal correlations other 
than flow are also discussed. 

Two particle correlations can also be used to make an 
event- by-event measurement of This method has two 
main advantages. The first advantage is the possibility 
to generate a mixed event background to calculate sta- 
tistical fluctuations. When single hits are mixed, the V2 
signal disappears whereas when pairs are mixed, the av- 
erage V2 signal is preserved. The second advantage is that 
there is only one event-by-event fit parameter, u^i since 
the reaction-plane dependence drops out when the dif- 
ference between the angles of particles is used. However, 
non-uniformities in acceptance are very hard to correct 
for with this approach, since the two-particle acceptance 
changes event-by-event with respect to the reaction-plane 
angle. 

Conceptually the simplest method is to measure V2 
event-by-event. In this method, absolute coordinates of 
the hits in the detector are used to measure V2 and the 



reaction plane. This approach also bears important dif- 
ficulties. Due to finite number fluctuations, the event- 
by-event V2 resolution is limited. As mentioned above, 
mixed events generated using single hits have no V2 signal 
and therefore cannot be used as a background reference. 
Furthermore, the resolution of the measurement changes 
with the true V2 value and the multiplicity in the event. 

In this paper, we will concentrate on the last approach. 
We will describe the method we have developed in or- 
der to address the difficulties outlined above for this ap- 
proach. 

In most fluctuation analyses, the statistical fluctua- 
tions are calculated using a mixed event background of 
certain event classes. In this analysis, events in the same 
event class correspond to events with the same reaction- 
plane angle. However the reaction-plane angle cannot 
be measured precisely, making it difficult to generate a 
mixed event background. Instead, MC simulations of the 
detector response will be used to account for statistical 
fluctuations. 

In a typical fluctuation measurement of a quantity x, 
the variance of x, can be decomposed into con- 

tributions from the statistical fluctuations, cr^stat: ^^'^ 
dynamical fluctuations, ct^ ^y^: 

2 2 I 2 fr-\ 

'^x,ohs — '''x.dyn + '^x,stat 

This equation holds if the average of the measurement, 
(0;°^*), gives the true average in the data, (x), and if the 
resolution of the measurement is independent of the true 
value. Neither of these conditions are satisfled in the 
event- by-event measurement of V2- Therefore, a more 
detailed study of the response function is required. 

We define K{v^^ tV2) as the distribution of the evcnt- 
by-event observed elliptic flow, Uj'"* , for events with con- 
stant input value of V2- If a set of events have an input 
V2 distribution given by /(W2), then the distribution of 
vf"^, g{vf'^), will be given by: 

g{vt)= / K{vf^,V2)f{v2)dv2 (6) 
Jo 

It is important to note that Eq. [5] holds in general and 
Eq.[5]can be derived from it in the special case described 
above. 

Thus, we separate the event-by-event elliptic flow fluc- 
tuation analysis into 3 tasks: 

• Finding givf"^) for a set of events, by an event-by- 
event measurement of . 

• Calculating the kernel, K{v2^'^ ,V2)j by studying the 
detector response. 

• Calculating the true V2 distribution, /(U2), by find- 
ing a solution to Eq. [6l 



4 



^ 




FIG. 3: The distribution of reconstructed hits in the 
PHOBOS muhiphcity array for events in a narrow vertex bin. 



B. Event-By-Event Measurement 

The PHOBOS detector employs silicon pad detectors 
to perform tracking, vertex detection and multiplicity 
measurements. Details of the setu p a nd the layout of 
the silicon sensors can be found in hj. The PHOBOS 
multiplicity array covers a large fraction of the full solid 
angle. At midrapidity, the vertex detector and the octag- 
onal multiplicity detector have different pad sizes and the 
acceptance in azimuth is not complete. Fig. [3] shows the 
distribution of reconstructed hits in the multiplicity ar- 
ray. The event-by-event measurement method has been 
developed to use all the available information from the 
multiplicity array to measure a single value, ^2^*, while 
allowing an efficient correction for the non-uniformities 
in the acceptance. 

The maximum likelihood method was applied for this 
purpose. We model the measured pseudorapidity depen- 
dence of V2 as: 



V2{v)^V2{v^0)-{l~\r^\/6) 



(7) 



This ansatz describes the main feature of the pseudora- 
pidity dependence of V2 over a range of centralities, shown 
in Fig. [31 We will denote V2{t] = 0) shortly as V2- We 
define the probability distribution function (PDF) of a 
particle to be emitted in the direction (77, (p) for an event 
with V2 and reaction plane angle (j)o- 

Pir,, q^\v2,M = p{ri)[l + 2v2iv) cos(2((/. - 0o))], (8) 

where p(r])d'i] is the probability for the particle to fall 
between 77 and i] + di] and [1 -I- 2v2{vi) cos(2(0 — ^o))]'^'/' 
is the probability of the particle to fall between (f) and 
(j) -\- d(f>. For a single event, the likelihood function of V2 
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FIG. 4: Elliptic flow as a function of pseudorapidity V2{rf) for 
charged hadrons for Au+Au collisions at ^/snn ~ 200 GeV 
for the three different centrality classes, 25-50%, 15-25% and 
3-15% central from top to bottomj^. The black lines show 
the fit given by Eq. [T] 



and do is defined as: 



(9) 



where the product is over all n hits in the detector. The 
likelihood function describes the probability of observing 
the hits in the event for the given values of the param- 
eters V2 and (1)0. Treating the hit positions as constants 
and the parameters V2 and c/io as variables allows choos- 
ing an estimate of the parameters, vf'" and (pQ^^ , which 
render the likelihood function as large as possible. When 
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comparing different values of V2 and 00 for a given event, 
p{ri) is constant and does not contribute to the calcula- 
tion. However, it is crucial that the PDF folded with the 
acceptance is normalized to the same value for different 
sets of parameters V2 and 0o- With these considerations, 
the PDF for hit positions in the detector is redefined. 
If the acceptance is given by A{ri,(f)), the normalization 
constant s is calculated in bins of rj as: 

X [l + 2t.2(l- h'l/e) cos(2(0-0o))]d0d?7' (10) 
Then the likelihood function becomes: 
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L{v2,(l>o) = Y[ 



1 



fj^ s(w2>o|'7») 
X [1 + 2t;2 (1 - cos(2(0, - (j)o))] (11) 

Instead of maximizing L{v2, (j^o), it is more convenient 
for technical reasons to maximize the auxiliary function 
K''^2,4>o), defined as: 

l{v2,(f>o) = V Inj — — -— 

X [l + 2«2(l- h.l/6) cos (2(0, -(/.o))]}, (12) 

where the sum is over the n reconstructed hits in the 
detector. Maximizing l{v2, 0o) as a function of V2 and 0o 
allows us to measure vf'" event-by-event. 

C. The Response Function 



To go from the measurement g(v2 ) to the true V2 
distribution, /(W2), one needs to determine the response 
of the measurement. The response function is the ker- 
nel in Eq. [S] Finite number fluctuations constitute the 
major part of the statistical fluctuations. Therefore, 
the kernel depends on the multiplicity in the detector: 
K — K , V2 , n) where n is the number of recon- 
structed hits. 

Monte Carlo simulations are used to 



measure 



K{vf'',V2,n). HIJING [18| is used to generate Au+Au 
events. The resulting particles in each event are redis- 
tributed in randomly with a probability distribution 
determined using Eq. [71 according to their 77 positions. 
These modified events are run through GEANT[l9j to 
simulate the PHOBOS detector response. 

When the event-by-event measurement is done on a 
set of MC events with a constant value of V2 and n, it 



FIG. 5: Distribution of the event-by-event measured quantity 
V2^" for MC events with a constant value of input V2 and a 
small range of hit multiplicity. Dashed line shows fit given by 
Eq.\M 



is observed that i^(w2^'*)|i,2,n is not Gaussian, but is well 
described by a Gaussian folded by a linear function: 



,,obs 



exp 



af 



262 



(13) 



as shown in Fig. [5l In the range of the measurement, 
the parameters a and b have a one-to-one correspon- 



dence with 



ohs\ 



and 



Therefore, the kernel can 



be found in MC simulations by measuring the distribu- 
tion of vf"' in bins of V2 and n and calculating {v^^) and 
(j^abs . We use 700,000 events over a collision vertex range 
of -10 cm< z <10 cm in this study, where z is the beam 
axis direction and the nominal vertex position is 2: = 0. 
The events are divided into 2 cm wide vertex bins. It is 
worth noting that the vertex resolution, which is better 
than 400 /im for a minimum-bias Au-|-Au event [2^, is 
much smaller than the vertex bins used in this analysis. 
{v2^^) and (J^obs are calculated in 70 bins in V2 and 8 bins 
in n. The results for one of the vertex bins are shown 
in Fig. [6l The statistical fluctuations in V2 and multi- 
plicity bins limit the precision of our knowledge of the 
kernel. Fitting smooth functions through the distribu- 
tions, (t'j'"') = {'"2'^){^2Tn) and a^obs = (j^obs{v2,n), re- 
duces the bin-by-bin statistical fluctuations and provides 
a simple parameterization of the kernel. The smooth flts 
to the results in Fig. [5] are plotted in Fig. [71 

The analysis is performed in bins of centrality. The 
paddle detectors which have pseudorapidity coverage of 
3.2 < 1 77 1 < 4.5 and cover 27r in the azimuthal direction 
are used in centrality determination. The details of the 
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FIG. 6; The mean (top) and the RMS (bottom) of the distri- 
bution of the event-by-event measured quantity V2"' for MC 
events in bins of input V2 and hit multiplicity. 



centrality determination can be found in d2.]. Assuming 
that the true V2 distribution, f{v2), is independent of the 
number of hits in the multiplicity array for a set of events 
in the same centrality class, it is possible to integrate out 
the multiplicity dependence according to: 



K{vf',V2)^ K{vf',V2,n)N{n)dn, (14) 



FIG. 7: Functions fitted through the distributions shown in 
Fig.H 



where N{n) is the distribution of number of hits for the 
set of events. This integration yields the appropriate ker- 
nel, plotted in Fig[9l for the given set of events with mul- 
tiplicity distribution N{n), shown in FiglHl 



D. Calculation of Dynamical Fluctuations 

As discussed in section lTl A| the last step of the analysis 
is to extract /(W2) from 5(^2'"*) ^-^d K{v2^^,V2)- There 
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FIG. 8: Distribution of the observed number of hits for a set FIG. 9: The kernel calculated for the set of events with mul- 
of events. tiplicity distribution shown in Fig. |H] according to Eq. 1141 



are many possible ways to address this problem. The 
approach depends on the question that we are trying to 
answer. In this case, we are interested in the mean and 
standard deviation of f{v2)- Therefore we assume an 
ansatz with two parameters: 



exp 



2< 



(15) 



f{v2) is only defined for V2 > 0. Therefore this ansatz is 
physically meaningful for the values, cry^/{v2) « 1. The 
validity of this ansatz, in particular for large Ov^l{v2) will 
be further studied. 

For given values of (1)2) and a^^, it is possible to take 
the integral in Eq. [6] and calculate the expected distribu- 
tion gexp{vf''\{v2),o^^). Comparing gexp{vf''\{v2),(Jv^) 
with the observation in data and minimizing the de- 
fined as: 



X 



values of 



E 



^bins 



(16) 



and (7^2 can be obtained. 



E. Verification of the Complete Analysis Chain 

The whole analysis procedure was tested using sets of 
MC events, selected from the events used to construct 
the kernel. Approximately 20,000 15-20% central events 
with a Gaussian distribution of V2 were used to construct 
each set. It should be noted that the functional form of 
the distribution of V2 in these sets is the same as our 
ansatz (Eg. llSp and the pseudorapidity dependence of V2 
is the same as in our PDF(Eq. [7]). Sample sets were se- 
lected with (U2) = 0.02, 0.03, 0.04 and 0.05. For each 
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FIG. 10: The true 112 distribution for a set of selected MC 
events with multiplicity distrubution shown in Fig.|8] Smooth 
line shows the parent Gaussian distribution. 
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FIG. 11: Distribution of the event-by-event measured quan- 
tity, V2''", for the set of events with the multiplicity distribu- 
tion in Fig. |8] and input «2 distribution in Fig. 1101 
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FIG. 12: Solid line shows the reconstructed «2 distribution for 
the set of events with input V2 distribution shown in dashed 
lines (also in Fig. I10|l and kernel shown in Fig. [O] 
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FIG. 13: Reconstructed (t„2 s-s a function of input cr^j for 
different sets of 15-20% central MC events. Plots are shown 
for input = 0.02 (top left), 0.03 (top right), 0.04 (bottom 
left) and 0.05 (bottom right). Statistical errors from combin- 
ing 10 different vertex bins are shown. 

value of (1)2), six different values of cr^^ (0.01, 0.015, 0.02, 
0.025, 0.03, 0.035) were chosen. Each set was divided 
in the 10 collision vertex bins, for which the kernels were 
constructed. The multiplicity distribution N(n), the ker- 



nel K{v2''^,V2), the input distribution /i„(w2) and the 
measured distribution g{v2^^) are plotted in Figures [51 [H 
[TUl and [TT] for a set of events from one vertex bin with 
{V2) = 0.05 and cr^j = 0.02. Minimizing the defined in 
Eq. [ini yields the estimate of the true {V2) and a^^ from 
the measurement. For this set, the values were found to 
be {V2) — 0.048 and a^^ = 0.023. The Gaussian distri- 
bution with these values is the estimate of the true V2 
distribution from the measurement. The distribution is 
plotted in Fig. [12] together with the true MC 1)2 ditribu- 
tion. 

Results obtained from different vertex bins were aver- 
aged. Fig. [Ol shows the combined results. It is seen that 
the method successfully reconstructs the input fluctua- 
tions down to {V2) ~ 0.03. 



III. CONCLUSION AND OUTLOOK 

The participant eccentricity model has increased the 
interest in measuring elliptic flow fluctuations. We have 
introduced a new analysis approach to perform this mea- 
surement. In this approach, all the available information 
from the PHOBOS multiplicity array is used to deter- 
mine V2 event-by-event. The response function of the 
event-by-event measurement, containing the contribu- 
tion of statistical fluctuations and detector effects is cal- 
culated using Monte Carlo simulations. Non-statistical 
fluctuations are extracted by unfolding the response func- 
tion from the distribution of the event-by-event measure- 
ment. Our approach has been tested on small sets of 
MC events, selected from the events used to calculate 
the response function. The input fluctuations are recon- 
structed successfully for {V2) > 0.03. 

The analysis can readily be applied to real data. The 
next step in the analysis will be to systematically study 
how the differences between HIJING MC and data in 
dN/drj, V2{ti) and azimuthal correlations other than flow 
influence the results. Different MC event generators and 
MC events with different input V2{rj) distributions will 
be used for this purpose. 
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